9 Shallow Flow
Shallow flow models are introduced to reflect a flowing free-surface continuum, whose characteristic horizontal scale \(L\) (in-flow direction) is much larger than its characteristic vertical extent \(H\) (flow-depth direction). Consequently, the shallowness parameter \(\epsilon := H/L\) is very small (\(\epsilon <<1\)). There are many examples for this type of flow regime, namely flow in layers of the atmosphere, oceanographic flow, river flow or gravity-induced surface flow such as landslides, rock flall or avalanches. You find these type of flows also in engineering contexts, for instance granular flow in process engineering.
9.1 Physical situation
We will derive the shallow flow model from general mass and momentum balance laws, in which the Cauchy stress tensor is kept in its general form form to allow for a large class of materials. In order to keep things simple, we will consider a two dimensional situation with the horizontal coordinate indicated by \(x\) and the vertical coordinate indicated by \(z\) (see Figure 9.1).
Note, that extension to a third dimension is straight forward. The system variables are then given by
\[ \begin{align*} \text{density:} \; & \rho(t, x,z)\\ \text{velocity:} \; & \mathbf v(t, x,z) = (u(t, x,z),w(t, x,z))^T \\ \text{Cauchy stress tensor:} \; & \mathbf \sigma(x,z) = \left( \begin{matrix} \sigma_{xx}(t, x,z) & \sigma_{xz}(t, x,z) \\ \sigma_{xz}(t, x,z) & \sigma_{zz}(t, x,z) \end{matrix} \right) \\ \text{basal topography:} \; & b(t,x)\\ \text{free surface:} \; & s(t,x) \end{align*} \]
Mass and momentum balance within the flowing material are given by
\[ \begin{align*} \nabla \cdot \mathbf v & = 0 && \tag*{🟦}\\ \partial_t (\rho \mathbf v) + \nabla \cdot (\rho \mathbf v \otimes \mathbf v) &= \nabla \cdot \sigma + \rho \mathbf g, && \tag*{🟧} \end{align*} \tag{9.1}\]
in which \(\mathbf g\) denotes the vector of gravitational acceleration.
We will consider flows that are driven by gravitational acceleration, such that it is convenient to assume the coordinate system to be rotated. In a situation, in which the basal topography is not planar, we assume the coordinate system to be inclined according to an average inclination angle \(\zeta\). Note, how in a rotated coordinate system that aligns with the inclined plane, the direction of the gravitational acceleration deviates from the \(z\) axis. We have
\[\mathbf g = [g_x, g_z]^T=g[\sin \zeta, -\cos \zeta]^T,\]
where \(g\) is the magnitude of the vector of gravitational acceleration. The regular Cartesian coordinate system is of course recovered by setting \(\zeta = 0\).
Componentwise, Equation 9.1 \(🟦\) and \(🟧\) read
\[ \begin{align*} \partial_x u + \partial_z w & = 0 &&\tag*{🟦}\\ \partial_t (\rho u) + \partial_x ( \rho u^2) + \partial_z (\rho u w) &= \partial_x \sigma_{xx} + \partial_z \sigma_{xz} + \rho g \sin \zeta && \tag*{🟧}\\ \partial_t (\rho w) + \partial_x (\rho u w) + \partial_z (\rho w^2) &= \partial_x \sigma_{xz} + \partial_z \sigma_{zz} - \rho g \cos \zeta && \tag*{🟨} \end{align*} \tag{9.2}\]
In this lecture, we will consider a physical regime, in which the density \(\rho\) can be assumed constant.
9.2 Dimensionless formulation
A scaling analysis helps to investigate the relative importance of the individual terms in this PDE. For shallow flows that are gravity induced, we choose the free fall time scale as the characterictic time scale
\[T := \sqrt{L/g}.\]
Considering \(T\), as well as characteristic vertical scale \(H\) and horizontal scale \(L\), we deduce two different velocity scales, the horizontal velocity scale
\[U := L/T = \sqrt{L g}\]
and the vertical velocity scale
\[W := H/T = \epsilon \sqrt{L g} (=\sqrt{ \epsilon H g}).\]
The characteristic scales of the system are hence given by \[ \begin{align*} x &= L \tilde x & \quad z &= H \tilde z & \quad t &= T \tilde t \\ u &= U \tilde u & \quad w &= W \tilde v & \quad \sigma &= \rho g H \tilde \sigma. \end{align*} \]
In the previous chapters we described in detail, how to translate a PDE system in its dimensionless form, such that details of the procedure can be omitted here. Note, that as before, the essential step is to transform the differentials, here \[ \begin{align*} \partial_x & = \partial_{\tilde x} /L & \quad \partial_z & = \partial_{\tilde z} /H & \quad \partial_t & = \partial_{\tilde t} /T. \end{align*} \]
Also as before, we will drop the \(\tilde \cdot\) for better readability. Unless otherwise stated, the variables are to be understood as their dimensionless counterparts.
The balance laws Equation 9.2 transform into
\[ \begin{align*} \partial_x u + \partial_z w & = 0 && \tag*{🟦} \\ \epsilon Fr^2 \left( \partial_t u + \partial_x u^2 + \partial_z (u w) \right) &= \epsilon \partial_x \sigma_{xx} + \partial_z \sigma_{xz} + \sin \zeta && \tag*{🟧} \\ \epsilon^2 Fr^2 \left(\partial_t w + \partial_x (u w) + \partial_z w^2 \right) & = \epsilon \partial_x \sigma_{xz} + \partial_z \sigma_{zz} - \cos \zeta, && \tag*{🟨} \end{align*} \tag{9.3}\]
in which we have introduced the dimensionless parameters \(\epsilon\) denoting the shallowness parameter and \(Fr\) denoting the Froude number:
\[ \begin{align*} \epsilon = H /L \quad \text{and} \quad Fr = u_r /\sqrt{g H} \end{align*} \]
9.3 Modeling evolving interfaces
Both the free-surface and - in its most general form - the basal topography are evolving with respect to space and time. This calls for a technique that allows us to capture these evolving interfaces. Here, we introduce a technique based on level set functions.
We assume the free surface and the basal topography to be given as functions \(s(\mathbf x,t)\) and \(b(\mathbf x,t)\) (see figure Figure 9.1). In order to derive the so-called kinematic boundary conditions we furthermore assume that material particles that are located at the surface stay at the surface, and that material particles that are located at the basal topography stay at the basal topography respectively.
Free-surface and basal topography can hence be expressed as the zero set of a level set functions \(F^{(i)}(t,x,z)\) in which \(i \in \{b,s\}\):
\[ \begin{align*} F^{(s)}(t,x,z) & = z - s(t,x) \tag*{🟦} \\ F^{(b)}(t,x,z) & = b(t,x) - z \tag*{🟧}. \end{align*} \]
Our above assumptions translates into a vanishing material derivative at the zero set
\[ \begin{align*} \left. \frac{D}{Dt} \, F_i(t, x,z) \right|_{F_i = 0} =0 \quad \text{for} \quad i \in \{b,s\} \end{align*} \]
Expansion of the first relation (\(i=s\)) yields
\[ \begin{align*} 0 & = \left. \left( \partial_t F_s(t, x,z) + \mathbf{u} \cdot F_s(t,x,z) \right) \right|_{F_s = 0} \\ & = \left. \left( - \partial_t s(t,x) - u(t,x,z) \partial_x s(t,x) + w(t,x,z) \right) \right|_{F_s = 0}, \end{align*} \]
into which we can substitute \(z=s(x,t)\):
\[ \begin{align*} \partial_t s(t,x) + u(t,x,s(x,t)) \partial_x s(t,x) &= w(t,x,s(x,t)). \end{align*} \]
This is known as the free-surface kinematic boundary condition:
\[ \begin{align*} \partial_t s + u^{(s)} \partial_x s &= w^{(s)}. \end{align*} \]
Note, that superscript \(^{(\ast)}\) indicates evaluation at \(z=\ast\). Similarly, we can derive a kinematic boundary condition for the basal topography:
\[ \begin{align*} \partial_t b + u^{(b)} \partial_x b &= w^{(b)}. \end{align*} \]
9.4 Modeling the interfaces’ stress states
9.4.1 Dynamic BC at the free surface
At the free surface, we assume a zero stress state, namely
\[\mathbf{\sigma} \cdot \mathbf n^{(s)} = \mathbf 0.\]
The outward pointing normal \(\mathbf n^{(s)}\) (see figure Figure 9.1) depends on the position and local orientation free-surface. It can be conveniently expressed in terms of of level set \(F^{(s)}\):
\[ \begin{align*} \mathbf n^{(s)} = \frac{\nabla F^{(s)}}{|\nabla F^{(s)}|} = \frac{1}{\sqrt{1+(\partial_x s)^2}} \left(\begin{matrix} -\partial_x s \\ 1 \end{matrix} \right). \end{align*} \]
In components, the zero stress state then reads
\[ \begin{align*} \sigma_{xz}^{(s)} - \sigma_{xx}^{(s)} \partial_x s & = 0 \tag*{🟦}\\ \sigma_{zz}^{(s)} - \sigma_{xz}^{(s)} \partial_x s & = 0, \tag*{🟧} \end{align*} \]
with its dimensionless counterpart reading
\[ \begin{align*} \sigma_{xz}^{(s)} - \epsilon \sigma_{xx}^{(s)} \partial_x s & = 0 \\ \sigma_{zz}^{(s)} - \epsilon \sigma_{xz}^{(s)} \partial_x s & = 0. \\ \end{align*} \]
9.4.2 Dynamic BC at the sliding surface
There are several options to close for the stress state at the sliding surface. Here, we will follow Hutter (2003) to assume a dry Coulomb friction relation, in which the friction induced shear stress \(\mathbf S^{(b)}\) is tangential to the bottom topography \(b(x,t)\) and opposing the surface tangential velocity. The magnitude of the dry Coulomb frictional resistance is proportional to the normal stress \(\mathbf N^{(b)}\) pointing in direction of \(\mathbf n^{(b)}\) (see figure Figure 9.2).
Dry Coulomb friction can hence be expressed as:
\[\mathbf{\sigma} \cdot \mathbf n^{(s)} = \mathbf S^{(b)} + \mathbf N^{(b)},\]
\[ \begin{align*} \sigma_{xx}^{(b)} \partial_x b - \sigma_{xz}^{(b)} & = -\left( \mu \frac{u^{(b)}}{| \mathbf v^{(b)}|} + \partial_x b \right) N \tag*{🟦}\\ \sigma_{xz}^{(b)} \partial_x b - \sigma_{zz}^{(b)} & = -\left( \mu \frac{w^{(b)}}{| \mathbf v^{(b)}|} - 1 \right) N \tag*{🟧} \end{align*} \]
in which \(N\) denotes the normal stress due to the material’s weight that depends on the average inclination \(\zeta\) and the local orientation of the basal topography \(b\). The scaled basal friction relation reads:
\[ \begin{align*} \epsilon \sigma_{xx}^{(b)} \partial_x b - \sigma_{xz}^{(b)} & = -\left( \mu \frac{u^{(b)}}{| \mathbf v^{(b)}|} + \epsilon \partial_x b \right) N \\ \epsilon \sigma_{xz}^{(b)} \partial_x b- \sigma_{zz}^{(b)} & = -\left( \epsilon \mu \frac{w^{(b)}}{| \mathbf v^{(b)} |} - 1 \right) N \end{align*} \]
The following relations results from subtracting the respective kinematic and dynamic boundary conditions and turn out to be useful later:
\[ \begin{align*} \partial_t (s(t,x)-b(t,x)) + \left[ u^{(z)} \partial_x z - w^{(z)} \right]_b^s &= 0 \tag*{🟦}\\ \left[ u^{(z)} \left(\partial_t z + u^{(z)} \partial_x z - w^{(z)} \right) \right]_b^s &= 0 \tag*{🟧} \\ \left[ \sigma_{xx}^{(z)} \partial_x z - \epsilon^{-1} \sigma_{xz}^{(z)} \right]_b^s &= \left( \epsilon^{-1} \mu \frac{u^{(b)}}{||\mathbf u^{(b)} ||} + \partial_x b \right) N \tag*{🟨} \end{align*} \tag{9.4}\]
9.5 Shallowness approximation and depth-integration
Assuming shallowness translate in assuming \(\epsilon <<1\). Neglecting terms that scale with \(\epsilon\) in the vertical momentum balance Equation 9.3 \(🟨\) yields
\[ \begin{align*} 0 &= \partial_z \sigma_{zz} - \cos \zeta, \end{align*} \]
mimicking a hydrostatic pressure distribution for isotropic material.
Ultimately, we are interested in the flow height, also referred to as flow depth and hence define the height of the flow
\[ \begin{align*} h(t,x):=s(t,x)-b(t,x). \end{align*} \]
We introduce the depth-weighted integral as
\[\bar \ast = (1/h) \int_b^s \cdot \,dz.\]
9.5.1 Depth-integrating the mass balance
We can now integrate the momentum balance Equation 9.3 \(🟦\) in \(z\)-direction:
\[ \begin{align*} & \int_b^s \partial_x u \,dz + \int_b^s \partial_z w \,dz = \partial_x ( h \bar u ) - \left[ u^{(z)} \partial_x z - w^{(z)} \right]_b^s = \\ & \partial_t h + \partial_x ( h \bar u ) = 0. \end{align*} \]
Note, that we made use of the Leibniz integral rule as well as Equation 9.4 \(🟦\), in which we identified \(h\).
9.5.2 Depth-integrating the momentum balance
Similarly we integrate the horizontal Equation 9.3 \(🟧\). Starting with the left hand side yields
\[ \begin{align*} & Fr^2 \int_b^s \left( \partial_t u + \partial_x u^2 + \partial_z (u w) \right) dz = \\ & Fr^2 \partial_t (h \bar u) + \partial_x (h \bar{u^2}), \end{align*} \]
in which we exploited Equation 9.4 \(🟧\).
Finally, we integrate the momentum equation’s right hand side
\[ \begin{align*} & \int_b^s \left( \partial_x \sigma_{xx} + \epsilon^{-1} \partial_z \sigma_{xz} + \epsilon^{-1} \sin \zeta \right) dz \\ & = \partial_x \int_b^s \sigma_{xx} dz - \left( \left[ \sigma_{xx}^{(z)} \partial_x z - \epsilon^{-1} \sigma_{xz}^{(z)} \right]_b^s \right) + \epsilon^{-1} h \sin \zeta \\ &= \partial_x ( h \bar \sigma_{xx} ) - \epsilon^{-1} \left( \mu \tfrac{u^{(b)}}{|| \mathbf u^{(b)} ||} + \epsilon \partial_x b \right) N_x + \epsilon^{-1} h \sin \zeta \end{align*} \]
and combine all into the depth-integrated mass and momentum balances
\[ \begin{align*} \partial_t h + \partial_x ( h \bar u ) &= 0 \\ Fr^2 \left( \partial_t (h \bar u) + \partial_x (h \bar{u^2} - Fr^{-2} h \bar \sigma_{xx} ) \right) &= - \epsilon^{-1} \left( \mu \tfrac{u^{(b)}}{|| \mathbf v^{(b)} ||} + \epsilon \partial_x b \right) N + \epsilon^{-1} h \sin \zeta. \end{align*} \]
9.5.3 Shape factor
It turns out that \(\bar{u^2}\) can be easily eliminated with th following trick: We write
\[ \bar{u^2} = \left. \bar u \right.^2 \underbrace{\frac{\bar{u^2}}{\left. \bar u \right.^2}}_{=\alpha} = \alpha \left. \bar u \right.^2 \]
and call \(\alpha\) the shape factor.
Rather than allowing \(\bar{u}\) and \(\bar{u^2}\) to be independent variables in our system, we describe \(\bar{u^2}\) according to its deviation from \(\bar{u}\). The shape factor \(\alpha\) can be calculated for standard parametrized velocity profiles (see exercise).
9.5.4 Constitutive closure
The depth-averaged stress \(\bar \sigma_{xx}\) is closed by considering that we are dealing with an isotropic material (such as water), in which the diagonal entries of the stress tensor reduce to negative pressure. The vertical momentum balance Equation 9.3 \(🟨\) then reads
\[ \begin{align*} 0 = \partial_z \sigma_{zz} - \cos \zeta = -\partial_z p - \cos \zeta \end{align*} \]
Subsequent depth-integration (again assuming a stress free state at the free surface) yields \[ \begin{align*} p(z) = -\int_b^z \cos \zeta \,dz + C = (h - (z-b))\cos \zeta \end{align*} \]
Integration furthermore yields \[ \begin{align*} h \bar \sigma_{xx} & = \int_b^s \sigma_{xx} dz = \int_b^s \sigma_{zz} dz = - \int_b^s p dz \\ & = - \int_b^s (h - (z-b))\cos \zeta \,dz = - \frac{h^2}{2} \cos \zeta. \end{align*} \]
9.6 Shallow flow model
Our findings so far can be summarized into a closed shallow flow model:
\[ \begin{align*} \partial_t h & + \partial_x ( h \bar u ) = 0\\ \partial_t (h \bar u) & + \partial_x \left(\alpha h \left. \bar u \right.^2 + Fr^{-2}\frac{h^2}{2} \cos \zeta \right) = \\ & - Fr^{-2} \epsilon^{-1} \left(\left( \mu \tfrac{u^{(b)}}{|| \mathbf u^{(b)} ||} + \epsilon \partial_x b \right) N + \epsilon^{-1} h \sin \zeta \right) \end{align*} \tag{9.5}\]
The system now comprises a set of two equations in one space dimension \(x\) as well as time \(t\), \(\alpha\) is the shape vector, \(Fr\) the Froude Number, \(h\) the height, \(\bar u\) the depth-averaged velocity, \(\mu\) the dry friction coefficient, \(\epsilon\) the shallowness coefficient, \(b(x)\) a function that denotes the basal topography, \(N\) the local normal force and \(\zeta\) the average inclination angle.
Many different flavours of shallow flow models exist that assume different friction relations, consitutive closures or even start from a two-phase system. A frictional relation that is for instance often found in geohazard engineering is a superposition of dry Coulomb friction with a so-called turbulent friction (Chezy friction) that scales with the square of the velocity.
In order to study the limiting regimes of the previously derived shallow flow model, we assume the flow model to be given in a simplified, generic form, in which
- we assume (almost) block flow, hence \(\alpha \equiv 1\), and consequently \(\partial_z u =0\)
- a very moderate average inclination, hence \(\zeta\) being small, and consequently \(\cos \zeta \approx 1\) and \(\sin \zeta \approx \zeta\)
- flow in downstream direction, hence \(u>0\)
- we assume a flowing material, for which a superposition of dry Coulomb friction (coefficient \(\mu\)) and turbulent friction (coefficient \(f\)) makes sense, e.g. snow.
- for the sake of convenience, we write \(u\) instead of \(\bar u\).
The shallow flow system Equation 9.5 reduces to
\[ \begin{align*} & \partial_t h + \partial_x ( h u ) = 0 \\ & Fr^{2} \epsilon \left( \partial_t (h u) + \partial_x (h u^2) \right) + \epsilon \partial (h^2/2) \\ & = \underbrace{\zeta h}_{\text{gravitational acceleration}} \; \underbrace{\underbrace{- \mu h}_{\text{dry Coulomb friction}} \; \underbrace{- f u^2}_{\text{turbulent Chezy friction}}}_{\text{friction relation}} \end{align*} \tag{9.6}\]
9.6.1 Long wave approximation
If we are interested in long range effects of a river it makes sense to choose a characteristic length scale according to the length of the river system, say \(L=10\)km, and a characteristic height according to the average depth of the river, say \(1\)m. With these scales we get a shallowness parameter which is \(\epsilon = H/L = 10^{-4}\). Neglecting terms that scale with \(\epsilon\) yields
\[ \begin{align*} \partial_t h + \partial_x ( h u ) &= 0\\ 0 &= \zeta h - \mu h- f u^2 \end{align*} \]
The momentum balance boils down to an algebraic equation in \(h\) and \(u\) (no derivatives involved), and can be solved explicitly for \(u\)
\[ \begin{align*} u = \sqrt{ \frac{\zeta - \mu}{f} h} \end{align*} \]
Here, we used our assumption that the material flows downstream (\(u>0\)). If we substitute this expression into the depth-averaged mass equation, we get the non-linear, scalar partial differential equation
\[ \begin{align*} \partial_t h + \sqrt{\frac{\zeta - \mu}{f}} \partial_x h^{3/2} = 0 \end{align*} \]
In order to get a feeling for this equation, we linearize it by applying the chain rule
\[ \begin{align*} \partial_t h + \frac{3}{2}\underbrace{\sqrt{\frac{\zeta - \mu}{f}} \sqrt{h}}_{= a(h)} \partial_x h = 0 \end{align*} \]
This hence constitutes a non-linear scalar advection equation with local advection speed \(a(h)\). Recall, that for a constant adevction speed \(a(h) \equiv a\), an initial wave profile \(h(t_0,x) = h_0(x)\) propagates at constant speed:
\[ \begin{align*} h(t,x) = h_0(x- a t). \end{align*} \]
As soon as the advection speed depends on \(h\) shocks can develop.
9.6.2 Short wave approximation
If to the contrary we consider a situation in which horizontal length scale is rather small in relation to the height, e.g. \(L=100\)m, vs \(H=1\)km, and hence \(\epsilon \approx 10\), we observe that \(\epsilon^{-1}\) tends to get small. Note that the depth-averaged equations are still ‘valid’ if we assume a vertically constant velocity profile. In that case Equation 9.6 reduces to
\[ \begin{align*} \partial_t h + \partial_x ( h u ) &= 0 &&\tag*{🟦}\\ Fr^{2} \left(\partial_t (h u) + \partial_x (h u^2) \right) + \partial_x (h^2/2) &= 0 && \tag*{🟧} \end{align*} \tag{9.7}\]
System Equation 9.7 is known as the shallow water equations.
The shallow water system is a famous PDE that has been successfully applied in many physical situations, for instance weather forecasteing. In the subsequent sections, we are interested in its properties.
9.7 Small purturbations or ‘a lake at rest’
We start by studying the onset of small perturbations based on the dimensional version of the shallow water system:
\[ \begin{align*} \partial_t h + \partial_x ( h u ) &= 0\\ \partial_t (h u) + \partial_x \left( h u^2 + g \,h^2/2 \right) &= 0 \end{align*} \tag{9.8}\]
Suppose, the shallow water system is at rest (\(u \approx 0\)) with uniform depth \(H\). This state trivially satisfies Equation 9.7. As known from the previous sections, we introduce a small perturbation as an asymptotic expansion: \[ \begin{align*} h(\mathbf x, t) &= \underbrace{h_0(\mathbf x, t)}_{=H \text{ "uniform depth"}} \quad + \quad \delta \; h_1(\mathbf x, t) + \quad ... && \tag*{🟦}\\ u(\mathbf x, t) &= \underbrace{u_0(\mathbf x, t)}_{=0 \text{ "at rest"}} \quad + \quad \delta \; u_1(\mathbf x, t) + \quad ... && \tag*{🟧} \end{align*} \tag{9.9}\]
in which \(\delta\) is the expansion coefficient. Substitution of the asymptotic expansions Equation 9.9 \(🟦\) and \(🟧\) into mass conservation of Equation 9.7 \(🟦\) yields \[ \begin{align*} \partial_t \left( H + \delta h_1 \right) + \partial_x \left( (H + \delta h_1) \, \delta \,u_1 \right) &= 0. \end{align*} \]
By expanding and exploiting the fact that \(H\) is constant reduces to
\[ \begin{align*} \delta \left( \partial_t h_1 + H \, \partial_x u_1 + \delta \left( h_1 u_1 \right) \right) &= 0. \end{align*} \]
All in all this can be combined into the following wave equation:
\[ \begin{align*} \partial^2_t h_1 = g H \partial^2_x h_1 \end{align*} \tag{9.10}\]
This is a wave equation that describes how small amplitude surface height perturbations evolve. Since we are interested in small perturbations it is natural to consider a solution of the form
\[ h_1 (\mathbf x,t) = A e^{\mathbf i (\kappa x - \omega t)} \tag{9.11}\]
Substituting Equation 9.11 into Equation 9.10 yields
\[ \begin{align*} \partial_t h_1 = - \mathbf i \omega h_1 \quad \text{and} \quad \partial^2_t h_1 = \omega^2 h_1 \\ \partial_x h_1 = \mathbf i k h_1 \quad \text{and} \quad \partial^2_x h_1 = \kappa^2 h_1 \end{align*} \]
Combining both yields
\[ \begin{align*} \omega^2 \pm g H \kappa^2 \end{align*} \]
The phase velocities are hence given by \(\omega/\kappa \pm \sqrt{g H}\). These are in particular constant, and hence we are facing a non-dispersive wave. In the context of the shallow water equations, \(\sqrt{g H}\) is called the celerity and determines the propagation speed of small amplitude perturbations.
9.8 Subcritical and supercritical flow
In this final section of the chapter on shallow flow, we want to study the different flow regimes of the non-linear shallow water system Equation 9.7. These can be characterized in terms of the system’s characteristic speeds, which are computed to be: \[ \begin{align*} \lambda_{1,2} = u \pm \sqrt{g h} = u \left( 1 \pm \frac{\sqrt{g h }}{u} \right) \end{align*} \]
The term \(\frac{\sqrt{g h}}{u}\) reminds us of the (inverse of the) previously introduced Froude number. Note however, that this time its components are not constant characteristic scales, but rather state variables that are varying with space and time. We therefore refer to \(\frac{u}{\sqrt{g h}}\) as the local Froude Number \(Fr_l(h,u)\). This allows us to discriminate two cases
\(Fr_l > 1 \quad \rightarrow \quad 1 + Fr_l^{-1} > 0 \; \text{and} \; 1 - Fr_l^{-1} > 0\),
hence both characteristic velocities point in flow direction. This regime is called supercritical or rapid.\(Fr_l < 1 \quad \rightarrow \quad 1 + Fr_l^{-1} > 0 \; \text{and} \; 1 - Fr_l^{-1} < 0\),
hence the characteristic velocities point in different directions. This regime is called subcritical or tranquil.
In either case, we observe that the shallow water system has real characteristic velocities (similar to the scalar advection equation). Recall that in these (hyperbolic) systems, information propagates at finite speeds, which has implications on the selection of an appropriate numerical solution scheme (FV or DG).